#!/usr/bin/python

import math
from math import log
from math import pow
from math import pi

from os import system

c=2.998e8# m/s
########################################################################
###### Parameters ######################################################
########################################################################
distance = 30000 # m

v0 = 0.7 # [light speed]
gamma = 1.79e4

########################################################################
########################################################################
########################################################################
names = list()
maxksi = 4
for ksi in range(0, maxksi+1):
	for field in ["E", "H"]:
		alpha = ksi/float(maxksi) * pi/2*0.9
		if field == "E":
			n = [-math.sin(alpha), 0, math.cos(alpha)]
			yunit = "E[V/m]"
		else:
			n = [0,1,0]
			yunit = "H[A/m]"
		obs = [distance*math.cos(alpha), distance*math.sin(alpha)]
		deg = round(alpha/pi*180)
		name = "alpha %i deg"%(deg)
		filename = "alpha_%ideg%s"%(deg, field)
		names.append((field, name, filename));
		
		xml = open("in.xml", "w")
		xml.write("<model>\n")
		xml.write("	<source>\n")
		t_start = 0
		xml.write("		<discharge t_start='%e'>\n"%(t_start))
		xml.write("			<current>\n")
		xml.write("				<shape type='expexp' alpha='1e4' beta='1e5' Q='3.87'/>\n")
		xml.write("				<delay type='exponential' v_c='%e' gamma='%e'/>\n"%(v0,gamma))
		#xml.write("			<decrease type='exponential' eta='0.7e-3'/>\n")
		xml.write("				<decrease type='none'/>\n")
		xml.write("			</current>\n")
		xml.write("<trajectory type='bezier3' \
		P0='[0,0,0]' P1='[0,0,5e3]' P2='[0,1.4e3,0.4e3]' P3='[0,1.4e3,5.4e3]' />")
		xml.write("		</discharge>\n")
		xml.write("	<mirror/>\n")
		xml.write("	</source>\n")
		xml.write("	<observer r='[%e,0,%e]' field='%s' n='[%e,%e,%e]'/>\n"%(obs[0],obs[1],field, n[0], n[1], n[2]))
		xml.write("</model>\n")
		xml.close()
		cfg=open("cfg.m", "w");
		cfg.write("ptxt = '%s'\n"%(name));
		cfg.write("yunit = '%s'"%(yunit))
		cfg.close()
		system("cat in.xml | ./qgr > %s.dat; cp in.xml %s.xml"%(filename, filename))


cfg=open("cfg.m", "w");
cfg.write("signals = {")
first = True
for row in names:
	if first:
		first = False
	else:
		cfg.write(";\n")
	cfg.write("'%s', '%s', '%s'"%(row[0], row[1], row[2]))
cfg.write("};\n")
cfg.close()
system("octave vis.m")
print "Done"


